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Abstract 

A  three  channel  polarimetric  deconvolution  algorithm  was  developed  to  mitigate  the  degrading  effects  of 
atmospheric  turbulence  in  astronomical  imagery.  Tests  were  executed  using  both  simulation  and  laboratory  data. 

The  resulting  efficacy  of  the  three  channel  algorithm  was  compared  to  a  recently  developed  two  channel  approach 
under  identical  conditions  ensuring  a  fair  comparison  amongst  both  algorithms.  Two  types  of  simulations  were 
performed.  The  first  was  a  binary  star  simulation  to  compare  resulting  resolutions  between  the  three  and  two  channel 
algorithms.  The  second  simulation  measured  how  effective  both  algorithms  could  deconvolve  a  blurred  satellite 
image.  The  simulation  environment  assumed  the  key  parameters  of  Fried’s  Seeing  parameter,  ,  and 

telescope  lens  diameters  of  and  .  The  simulation  results  showed  that  the  three  channel  algorithm 

always  reconstructed  the  true  image  as  good  as  or  better  than  the  two  channel  approach,  while  the  total  squared  error 
was  always  significantly  better  for  the  three  channel  algorithm.  The  next  step  is  comparing  the  two  algorithms  in  the 
laboratory  environment.  However,  the  laboratory  imagery  was  not  actually  blurred  by  atmospheric  turbulence,  but 
instead  camera  defocusing  was  used  to  simulate  the  blurring  that  would  be  caused  by  atmospheric  turbulence.  The 
results  show  that  the  three  channel  significantly  outperforms  the  two  channel  in  a  visual  reconstruction  of  the  true 


image. 
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1.  Introduction 

1.1  Motivation  for  Research  and  Applications 

Polarimetric  imaging  in  astronomy  offers  one  distinct  advantage  when  compared  to  standard  imaging 
techniques.  This  advantage  arises  when  imaging  man-made  objects,  because  of  the  tendency  for  these  objects  to 
reflect  polarized  light  more  than  its  surroundings.  The  polarized  data  that  is  collected  from  objects  in  space  can  be 
very  beneficial  in  image  analysis,  and  is  proven  in  this  research.  However,  polarimetric  imaging  performance  is 
limited  by  the  blurring  caused  by  the  atmosphere.  Therefore,  it  is  very  beneficial  to  combine  polarimetry  with  an 
image  processing  technique  known  as  deconvolution  to  deblur  images  degraded  by  the  atmospheric  conditions. 

One  of  the  main  goals  of  this  research  is  to  develop  a  three  channel  polarimetric  deconvolution  algorithm 
with  superior  resolving  capability  to  that  of  the  prior  researched  two  channel  approach  [6].  A  two  point  resolution 
metric  is  used  in  this  research  to  determine  which  imaging  technique  is  superior.  The  technique  is  similar  to  that 
used  by  Lord  Rayleigh  in  that  resolution  is  determined  by  the  ability  to  resolve  two  point  sources  [3].  A  majority  of 
simulations  and  critiques  to  the  algorithm  proposed  in  this  thesis  were  centered  on  the  idea  of  resolving  two  point 
sources.  The  efficacy  of  the  deconvolution  algorithm  to  resolve  two  point  sources  is  directly  related  to  how  well  it 
will  deconvolve  a  more  complex  image,  such  as  a  satellite  image. 

Another  goal  of  this  research  is  to  not  only  resolve  sharper  images,  but  also  to  collect  and  estimate 
polarization  data  from  an  object  in  observation.  Polarization  data  can  be  used  specifically  to  assist  in  identifying  a 
target,  or  more  generally  to  identifying  materials  [5]. The  proposed  algorithm  is  capable  of  estimating  a  parameter 
defined  as  the  polarization  attenuation  factor.  This  factor  is  mathematically  dependent  on  the  angle  of  polarization  of 
the  light  reflected  by  the  material  in  the  observation.  The  angle  of  polarization  data  is  the  potentially  useful 
information  in  identification  of  materials  or  targets. 
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The  capability  to  analyze  sharper  astronomical  images  and  corresponding  polarization  data  can  prove  to  be 
very  useful  for  the  United  States  Air  Force.  As  the  DoD’s  perception  of  war  in  space  evolves  it  is  necessary  that  we 
always  maintain  space  superiority.  Cutting  edge  surveillance,  intelligence,  and  reconnaissance  are  imperative  to 
attaining  space  superiority,  and  that  is  where  this  research  potentially  contributes.  For  example,  someday  it  may  be 
very  useful  to  closely  monitor  our  satellites  and  enemy  satellites.  Satellites  in  extreme  distant  orbits  (Geostationary) 
do  not  return  much  light  and  provide  very  dim  and  low  resolution  images.  If  an  enemy  satellite  was  in  close 
proximity  to  our  own  satellite,  we  may  be  interested  to  know  that.  Unfortunately,  an  image  from  current  telescopes 
may  not  show  a  distinction  between  two  satellites  in  relatively  close  proximity,  but  more  likely  will  blur  them 
together  into  one  mass  of  light.  However,  if  polarimetry  in  conjunction  with  deconvolution  is  used  it  may  become 
more  apparent  that  there  are  two  satellites  instead  of  one.  The  deconvolution  provides  for  more  clarity  in  the  image 
by  deblurring,  while  the  polarization  information  can  help  in  determining  if  the  object  in  observation  is  man-made. 
This  type  of  capability  is  just  one  example  of  how  the  USAF  can  retain  space  superiority. 

1.2  Problem  Statement 

One  major  problem  with  imaging  objects  in  space  is  the  negative  effects  the  atmosphere  introduces  into  the 
imaging  process.  Specifically,  atmospheric  turbulence  causes  aberrations  in  images  when  viewed  from  the  surface  of 
earth  [12:69].  Figure  1.1  demonstrates  the  problem  presented  by  atmospheric  turbulence.  This  scenario  shows  a  star 
that  is  a  great  distance  away  from  the  earth.  The  light  from  the  star  will  naturally  radiate  in  a  spherical  manner,  but 
as  it  travels  a  distance  of  “near  infinity”  it  becomes  a  plane  wave.  As  the  wave  enters  the  earth’s  atmosphere  it  starts 
out  flat,  but  as  it  travels  through  the  phase  front  changes  and  becomes  degraded. 

Atmospheric  turbulence  causes  time  delays  in  the  light  wave  as  it  passes  through  the  atmosphere,  but  the 
delay  is  not  uniform  across  the  entire  plane  wave.  The  tilt  illustrates  how  the  right  portion  of  the  light  passes  through 
after  the  left  side.  This  delay  will  inevitably  have  an  effect  on  the  image  formed  at  the  detector  of  camera,  which  is 
on  the  other  side  of  the  lens  of  a  telescope.  The  lens  of  the  telescope  will  cause  all  incoming  light  outside  to 
converge  as  shown  towards  the  detector.  An  ideal  scenario  would  be  that  the  light  comes  in  untilted  and  forms  an 
image  centered  on  the  detector,  but  in  reality  this  is  not  the  case.  An  input  of  tilted  light  on  a  lens  will  result  in  a 
shifted  image  on  the  corresponding  detector,  based  on  the  shift  theorem  of  Fourier  transforms  [3:8].  After  a  “long” 
exposure  time,  there  will  not  be  just  one  shift  of  the  point  source,  but  an  accumulation  of  many  shifts,  which  appear 
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as  blur  in  the  image  There  are  other  aberrations  in  the  imagery  due  to  atmospheric  effects,  but  tilted  light  is  the 
dominant  source  of  blurring  in  long  exposure  astronomical  imagery  [12:69].  Image  A  in  Figure  1.1  shows  a 
simulated  example  of  an  ideal  image  of  two  point  sources  (top  of  figure),  and  the  resulting  blurry  measured  image  at 
the  detector  (Image  B).  The  bottom  image  in  the  figure  will  be  referred  to  from  now  on  as  measured  data  d(x,y). 


Figure  1.1:  Left-Instantaneous  example  of  the  consequences  of  atmospheric  turbulence  on  a  traditional  non- 
polarimetric  imaging  system.  Image  A-  Example  of  two  binary  star  true  intensity.  Image  B-  Example  of  data 
measured  on  camera  detector  known  as  measured  data 


1.3  Scope  of  Research  Effort 


The  goal  of  this  research  is  to  provide  an  effective  three  channel  polarimetric  deconvolution  algorithm  for 
reconstructing  images  which  have  been  degraded  by  atmospheric  turbulence  and  noise.  This  thesis  explains  in 
entirety  the  process  used  for  deblurring  and  de-noising  images.  Results  are  shown  for  both  simulation  and  laboratory 
data. 
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1.4  Chapter  Summaries 

Chapter  1:  Introduction 

This  chapter  introduces  the  topic  of  three  channel  polarimetric  imaging,  by  first  explaining  the  motivation 
and  applications  of  it.  It  then  gives  a  concise  problem  statement,  followed  by  the  scope  of  this  actual  research  effort. 

Chapter  2:  Background 

Chapter  2  provides  a  chronological  background  of  all  the  research  performed  in  the  last  four  decades, 
which  have  significantly  contributed  to  the  research  effort  described  in  this  thesis. 

Chapter  3:  Simulation  Description 

A  description  of  the  simulation  is  the  most  important  part  of  this  thesis.  Chapter  3  starts  by  introducing  the 
convolution  model,  which  addresses  how  an  ideal  image  is  degraded.  The  convolution  model  is  broken  up  into  three 
subsections:  an  explanation  of  optical  transfer  functions,  the  noise  model,  and  a  section  to  explain  how  an  image  is 
actually  degraded  via  a  process  explanation.  The  second  subsection  to  Chapter  3  is  the  deconvolution  model,  which 
explains  the  process  for  inverting  the  convolution  model.  It  starts  by  introducing  the  process  of  deconvolution,  and 
then  transitions  into  a  more  technical  discussion  of  the  equations  and  derivation  necessary  to  have  a  functional 
deconvolution  algorithm.  It  finishes  with  a  brief  discussion  on  dual  point  source  pixel  spacing  for  the  ideal  image  in 
simulation. 

Chapter  4:  Simulation  Results 

The  simulation  results  section  is  broken  up  into  two  major  parts.  The  first  section  addresses  the  simulation 
results  for  dual  point  sources,  and  the  second  discusses  the  results  for  simulations  executed  with  an  ideal  satellite 
image.  It  is  shown  in  all  simulation  cases  that  the  three  channel  algorithm  outperforms  or  performs  the  same  as  the 
two  channel  in  aesthetics,  while  the  three  channel  always  outperforms  in  total  squared  error  calculations. 

Chapter  5:  Laboratory  Imagery  Results: 
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Laboratory  imagery  is  the  final  step  taken  in  this  thesis  to  compare  both  algorithms.  This  section  explains 
the  laboratory  set  up,  image  degradation  model,  and  analysis  of  results.  The  results  show  that  the  three  channel 
algorithm  outperforms  the  two  channel  algorithm. 

Chapter  6:  Limitations 

The  limitations  chapter  explains  potential  limitations  of  the  three  channel  algorithm  either  in  simulation 
scenarios  or  laboratory  data  scenarios.  The  four  limitations  explained  are  the  effect  of  lens  diameter  to  seeing 
parameter  (  )  ratio,  the  absence  of  blind  deconvolution,  lack  of  algorithmic  autonomy,  and  diversity  of  multi¬ 
channel  data  due  to  polarization. 
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2.  Background 


William  Richardson  began  the  modem  era  of  deconvolution  in  image  processing  when  he  published  an 
algorithm  that  was  dependent  on  Bayes’s  Theorem.  He  describes  in  his  journal  article  Bayesian-Based  Iterative 
Method  of  Image  Restoration  [11]  how  an  image  degraded  by  a  point  spread  function  can  be  restored  to  a  near 
perfect  original  image  all  through  use  of  Bayes’  theorem.  He  built  the  framework  for  many  future  researchers  when 
he  derived  an  iterative  approach  of  applying  Bayes’s  theorem  to  estimate  the  tme  values  of  an  image  given  the 
degraded  image  and  a  known  point  spread  function.  In  1974  Lucy  also  devised  a  clever  use  of  Bayes’s  theorem  in 
application  to  estimation  theory  and  he  took  it  one  step  further  by  including  the  concepts  of  maximum  likelihood  in 
his  estimation  process  [10].  Lucy’s  work  and  concepts  were  applied  to  astronomy,  but  not  astronomical  imaging, 
which  is  the  topic  of  this  research  effort. 

In  1976  Dempster,  Laird  and  Rubin  wrote  a  journal  article  for  the  Royal  Statistical  Society  called. 
Maximum  Likelihood  from  Incomplete  Data  via  the  EM  Algorithm  [1],  which  had  a  tremendous  impact  in  image 
processing  theory.  In  it  they  describe  a  process  called  the  Expectation  Maximization  (EM)  process.  The  paper 
presents  an  iterative  approach  to  maximum  likelihood  estimation  with  the  observed  data  being  referred  to  as 
incomplete  data.  The  EM  process  is  dependent  on  two  types  of  data,  the  first  being  incomplete  data  and  the  other 
complete  data.  The  complete  data  is  an  abstract  idea  and  its  relevance  and  application  to  this  research  is  explained 
later  in  3.2.2.  The  process  was  termed  the  EM  process  because  each  iteration  has  two  key  steps.  The  first  is  the 
expectation  step,  and  the  second  is  the  maximization  step  [1].  These  steps  are  also  more  clearly  defined  with  respect 
to  this  research  in  3.2.2. 

In  1982  Shepp  and  Vardi  wrote  an  article  on  emission  tomography,  where  they  used  concepts  of  maximum 
likelihood  estimation  and  EM  to  reconstruct  images.  This  document  is  a  major  source  of  interest  for  this  research 
because  of  the  developments  these  two  men  made  in  application  of  the  EM  process  [14].  Specifically,  their  results 
for  how  to  define  the  conditional  expectations  were  used  later  by  Timothy  Schulz  in  his  research  [13].  The 
mathematical  model  developed  throughout  this  thesis  is  most  closely  related  to  the  work  described  by  Schulz  in  his 
1993  Joint  Optical  Society  of  America  article  [13],  Section  3.3.2  details  the  mathematical  model  used  in  this 
research  and  a  detailed  discussion  of  the  similarities  and  differences  between  Schulz’s  model  and  our  own. 
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About  thirteen  years  later  in  2006  Major  Dave  Strong  continued  researching  in  the  application  of  the  EM 
process  to  astronomical  imaging.  He  advanced  Schulz’s  research  one  step  further  and  developed  a  two  channel 
polarimeter,  which  used  one  channel  of  unpolarized  light  and  one  channel  of  linearly  polarized  light.  His  solution 
was  also  restricted  to  short  exposure  imagery  as  was  Schulz’s  [15]. 

Soon  after  Strong  finished  his  dissertation,  Daniel  LeMaster  researched  in  the  field  of  polarimetric  imaging. 
In  2007  LeMaster  researched  a  similar  three  channel  approach  to  what  has  been  developed  in  this  research,  but  there 
are  some  critical  differences  [8],  First  and  most  importantly,  he  was  only  researching  deconvolution  algorithms  for 
short  exposure  imagery,  where  this  effort  only  considered  long  exposure  imagery.  LeMaster  also  used  different 
polarization  analyzer  settings  (0°,  60°,  -60°)  than  used  in  the  three  channel  approach  described  throughout  this  thesis 
(0°,  45°,  90°).  The  difference  in  polarization  analyzer  settings  determines  a  different  derivation  resulting  in  a 
different  deconvolution  algorithm.  He  also  directly  estimated  the  angle  of  polarization  of  the  light  in  his 
deconvolution  algorithm,  where  he  discovered  that  the  estimates  of  the  angles  of  polarization  are  very  susceptible  to 
noise  and  time  consuming  to  obtain.  The  deconvolution  algorithm  developed  in  this  research  does  not  estimate  the 
angle  of  polarization,  but  it  does  estimate  the  polarization  attenuation  factors,  which  are  dependent  on  the  angle  of 
polarization.  Estimating  this  parameter  may  prove  to  be  less  susceptible  to  noise  and  require  fewer  mathematical 
steps.  This  concept  is  further  explained  in  3.1.3  and  3.2.2. 

LeMaster  was  able  to  effectively  restore  images  using  the  developed  three  channel  deconvolution 
algorithm,  but  he  was  never  able  to  effectively  estimate  the  angle  of  polarization.  LeMaster’s  dissertation  [8]  also 
lacked  one  important  thing  and  that  was  a  comparison  of  a  three  channel  versus  two  channel  algorithm.  Conversely, 
this  research  effort  is  focused  on  a  comparison  of  three  channel  versus  two  channel  algorithm.  The  two  channel 
algorithm  compared  to  was  developed  by  Captain  Steve  James  in  2008  [6].  In  his  research  he  developed  a  two 
channel  polarimeter  deconvolution  algorithm  capable  of  restoring  images  corrupted  by  noise  and  blurring.  His 
algorithmic  solution  was  determined  through  use  of  a  Bayesian  estimation  method  called  maximum  a  posteriori 
(MAP)  estimation  [7:  350].  He  also  utilized  the  concepts  of  complete  and  incomplete  data,  which  are  defined  in  the 
EM  process  [1],  [7:  187],  It  was  through  the  use  of  MAP  estimation  and  the  EM  algorithm  that  he  was  able  to 
develop  a  deconvolution  algorithm  capable  of  restoring  images  to  a  less  noisy  and  blurred  state.  The  three  channel 
approach  described  in  this  thesis  closely  mimics  his  strategy  of  using  MAP  estimation  in  conjunction  with  the  EM 
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process,  but  requires  a  different  mathematical  derivation  to  arrive  at  the  update  equations  inside  the  deconvolution 


algorithm. 
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Dempster,  Laird  and  Rubin 

Maximum  Likelihood  from  Incomplete  Data  via  the  EM  Algorithm 

1982 

Shepp  and  Vardi 

Maximum  Likelihood  Reconstruction  for  Emission  Tomography 

1993 

Schulz 

Multiframe  Blind  Deconvolution  of  Astronomical  Images 

2007 

Strong 

Polarimeter  Blind  Deconvolution  Using  Image  Diversity 
2008 
Lemaster 

Statistical  Processing  Methods  for  Polarimetric  Imaging 

2009 

James 

Blind  Deconvolution  Through  Polarization  Diversity  of  Long  Exposure  Imagery 

2011 

Engelson 

3  Channel  Polarimetric  Based  Data  Deconvolution 


Figure  2.1:  Literature  Review  Timeline:  chronologically  displays  contributions  used  to  develop  and  execute  3 
channel  polarimetric  deconvolution  of  astronomical  images. 
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3.  Simulation  Description 


3.1  Convolution  Model 

This  section  addresses  the  transformation  process  of  ideal  true  intensity  to  the  blurred  and  noisy  measured 
data.  First,  the  idea  of  true  intensity  must  be  understood.  One  example  of  true  intensity  can  be  seen  in  Image  A  of 
Figure3.1,  where  there  is  an  image  of  two  binary  stars.  It  is  called  true  intensity  because  it  is  a  perfect  and 
unaberrated  image.  It  is  also  the  benchmark  data  set  used  to  compare  results  to  later  on. 

Figure  3.1  demonstrates  the  process  by  which  the  true  intensity  is  processed  by  a  camera.  It  also  provides  a 
block  diagram  of  the  optical  system,  where  the  input  is  true  intensity  convolved  with  the  system  point  spread 
function  (PSF)  .  Due  to  the  properties  of  linearity  and  shift  invariance  the  system  can  be 

characterized  by  the  PSF,  as  an  electrical  system  can  be  characterized  by  its  impulse  response  [3: 19].  The  system 
PSF  is  a  combination  of  the  long  exposure  PSF,  ,  convolved  with  optics  based  PSF,  (refer  to  Equation  (3.4)). 
The  measured  data  is  the  result  of  a  convolution  of  true  intensity  with  plus  noise.  The  blurry  but  non- 

noisy  image  is  related  to  the  true  intensity,  via  equation: 


(3.1) 


It  is  space  invariant  because  only  depends  on  the  distances  between  points  and 

where  is  the  detector  plane  coordinates  and  is  the  observation  plane  coordinates  [3:18-21]. 
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Figure  3.1:  A  continuation  of  Figure  1  with  the  deconvolution  model  boxed  off  to  the  left.  Image  A-  True  intensity. 
Image  B-  True  intensity  blurred  with  long  exposure  transfer  function  .  Image  C-  True  intensity  blurred  with 
optics  based  OTF  .  Image  D-  True  intensity  with  added  Poisson  modeled  noise.  Image  E-  Data  measured  at 
camera  detector 


3.1.1  Optical  Transfer  Function  Model 

There  are  two  sources  of  blurring  specifically  considered  in  this  research  in  the  degradation  process  of  the 
data.  Image  B  in  Figure  3.1  demonstrates  the  first  source  of  blurring  is  caused  by  atmospheric  turbulence.  The 
blurring  caused  by  atmospheric  turbulence  is  modeled  with  the  long  exposure  optical  transfer  function  (LEOTF), 
which  is  the  Fourier  transform  of  the  PSF  .  It  is  the  cause  for  the  most  severe  degradation  to  the  true  intensity. 
This  OTF  is  only  used  in  long  exposure  imagery  (  ).  If  short  exposure  imagery  were  relevant  to  this 

research  as  it  was  to  some  of  the  previous  efforts  discussed  in  section  2,  then  the  OTF  would  be  modeled  with  a 
different  equation  than  the  LEOTF.  The  formula  for  the  LEOTF  is  given  by  Equation  (3.2)  [2:402,439].  Notice  that 
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the  resulting  magnitude  of  is  solely  dependent  on  r0  (Fried’s  Seeing  [2:429-431]).  Fried’s  Seeing  Parameter 
was  developed  to  assist  in  modeling  the  atmospheric  effects  in  imagery.  A  key  ratio  is  in  the  exponent  of  Equation 

(3.2),  where  the  ratio  describes  the  effect  of  the  atmosphere  in  the  image.  In  this  equation  is  the  mean 

wavelength  of  the  light,  is  the  spatial  frequency  and  is  the  effective  focal  length  of  the  optical  system.  If  the  ratio 
is  greater  than  1,  then  the  atmosphere  will  start  to  introduce  blurring  into  the  image.  If  it  is  less  than  1,  then  the 
resolution  of  the  image  is  unaffected  by  the  atmosphere.  [2:429-431].  Variables  and  represent  the  spatial 
frequencies  in  the  x,  y  plane  (camera  detector  plane),  where  the  OTF  is  only  relevant  in  the  spatial  frequency 
domain,  and 


(3.2) 

The  least  significant  contribution  in  the  blurring  of  the  true  intensity  is  due  to  the  optics  based  OTF  (  ) 

given  in  Equation  (3.3)  [2:364],  is  the  frequency  domain  representation  of  the  system  impulse  response  or  PSF, 
represents  the  pupil  function,  where  are  spatial  coordinates  in  units  of  meters.  In  this  research 
the  pupil  function  is  defined  at  the  lens  of  the  telescope,  where  the  pupil  function  directly  determines  the  resulting 
optical  transfer  function  [3:141].  The  pupil  function  accounts  for  the  effects  caused  by  an  aperture  in  a  optical 
system  [3:107].  Equation  (3.3)  is  only  true  for  the  case  of  an  incoherent  light  source.  The  blurring  caused  by 

accounts  for  diffraction  and  other  aberrations  in  the  measured  data.  In  this  research  diffraction  is  the  only  cause 
for  any  blur  to  the  true  intensity  with  respect  to  .  A  diffraction  pattern  can  be  seen  in  Image  C  of  Figure  3.1.  This 
example  shows  how  miniscule  the  diffraction  based  blurring  is  compared  to  the  LEOTF  blurring.  The  product  of 
both  OTFs  gives  the  OTF  for  the  system  shown  in  Equation  (3.3).  is  the  composite  of  both  OTFs  and  is 

used  for  degrading  all  true  intensity. 


(3.3) 


(3.4) 
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3.1.2 


Noise  Model 


There  are  two  types  of  noise  to  consider  in  this  research  and  those  are  readout  noise  and  photon  counting 
noise.  Both  types  of  noise  are  present  in  photon  counting  cameras.  Readout  noise,  also  referred  to  as  thermal  noise  is 
caused  by  imperfections  in  the  analog  to  digital  conversion  process  in  a  camera,  and  the  injection  of  currents  from 
other  electronics  in  the  circuitry.  Readout  noise  variance  can  be  measured  by  blocking  light  from  entering  the 
camera.  The  second  type  of  noise  considered  is  photon  count  noise,  and  is  assumed  to  be  the  dominant  source  of 
noise  in  this  research.  It  is  safe  to  assume  this  if  an  avalanche  photo-diode  (APD)  is  used.  The  thermal  noise 
becomes  negligible  due  to  the  high  amplification  effects  of  this  device.  In  this  case  the  signal  generated  by  each 
photon  is  amplified  to  a  level  above  the  readout  noise  level.  Photon  count  noise  is  caused  by  incoherent  light  (i.e. 
sunlight)  and  can  be  modeled  with  Poisson  statistics  [12:16],  [2:485]  .The  random  arrival  of  photons  counted  by  a 
camera  detector  is  the  consequence  of  photon  count  noise.  Figure  3.2  shows  an  integration  time  window  where  there 
are  eight  expected  photons  arrivals  along  with  the  actual  amount  of  nine  arrivals.  The  irregular  time  of  arrival  of  each 
photon  contributes  to  the  photon  counting  noise. 

Photon  counting  noise  is  best  described  by  use  of  the  Poisson  PMF  [2:485].  Equations  (3.5)  and  (3.6)  show 
the  general  probability  mass  function  for  a  Poisson  random  variable,  and  a  Poisson  PMF  relevant  in  notation  to  this 
research  [9:508].  In  the  general  Poisson  PMF  equation  is  the  random  variable  representing  the  actual  number  of 
occurrences,  is  the  mean  number  of  occurrences.  is  the  probability  that  the  number  of  event  occurrences 

is  equal  to  .  Similarly,  Equation  (3.6)  presents  a  two-dimensional  case  where  the  mean  number  of  occurrences 
is  ,  also  known  as  the  intensity,  and  the  actual  number  of  occurrences  (photon  counts  at  each  pixel)  is 


(3.5) 


(3.6) 
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Equation  (3.7)  shows  how  the  SNR  is  calculated  for  all  scenarios.  It  also  shows  how  one  of  the  two 
variables  in  the  denominator  can  be  the  dominant  source  of  noise,  if  it  is  much  greater  than  the  other.  The  gain  factor 
Gapd  of  the  APD  is  shown  in  bold  and  only  relevant  if  an  APD  is  used.  represents  the  thermal  noise  variance  in 
electrons.  is  the  signal  strength  in  photons  with  photon  counting  noise  included  in  the  signal  [12:23-24]. 

Assume  that  is  negligible  compared  to  the  noise  included  in  multiplied  by  GApd-  With  this  assumption 

the  SNR  ends  up  being  the  variance  of  the  signal  divided  by  the  standard  deviation  of  the  signal,  due  to  the 
properties  of  Poisson  random  variables. 


(3.7) 


Figure  3.2:  Random  Photon  Arrivals:  There  are  expected  to  be  8  photon  arrivals  in  the  integration  time  window, 
but  due  to  randomosity  of  photon  arrivals  there  actually  9  photon  arrivals,  and  not  at  expected  times. 

3.1.3  Simulating  Image  Degradation 


This  section  describes  the  simulation  process  for  degrading  the  ideal  true  intensity  to  what  is  referred  to  as 
the  measured  data  .  An  example  of  measured  data  can  be  seen  in  image  E  of  Figure  3.1.  Figure  3.3  displays 

an  example  diagram  of  a  potential  three  channel  polarimeter  setup.  Notice  that  each  polarization  analyzer  is  set  in  a 
particular  direction  and  that  there  are  three  detectors  for  collecting  three  channels  worth  of  “measured  data.”  The 
three  polarization  analyzers  combined  with  three  channels  of  measured  data  deem  the  necessary  title  3  Channel 
Polarimeter. 
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Polarization 

Analyzers 


Figure  3.3:  An  example  of  the  setup  for  a  three  channel  polarimeter  relevant  to  this  research  effort.  There  is  an 
assumed  incoherent  light  input  to  the  optical  system.  There  must  be  three  polarization  analyzers  set  as  specified 
above,  with  three  detectors  to  collect  data  from  each  analyzer. 

Figure  3.4  below  shows  the  logical  flow  of  the  data  simulation  process.  The  first  three  rhombuses  represent 
data  that  is  either  created  in  the  simulation  program  or  uploaded  for  use  in  the  simulation  program.  The  true  degree 
of  polarization  represents  the  degree  to  which  light  is  polarized  in  the  scene.  Only  certain  parts  of  the  object  in 
observation  will  reflect  polarized  light.  As  the  polarized  light  passes  through  the  polarization  analyzers  a  fraction  of 
the  light  will  be  lost,  unless  it  is  in  perfect  alignment  with  one  of  the  polarization  analyzers. 

The  data  represents  the  angle  of  polarization  of  each  element  in  the  polarization  data.  More 

importantly,  is  used  to  calculate  three  values,  which  corresponds  to  the  three  polarization 

attenuation  factors.  The  variables  ,  ,  and  indicate  the  relative  alignment  of  the  polarized 

light  with  each  channels’  polarization  analyzer. 

The  next  three  rectangles  represent  calculations  for  ,  ,  ,  LEOTF,  and  optics  based 

OTF.  Rectangle  four  initiates  the  process  for  degrading  the  true  intensity.  There  are  four  sets  of  data  that  result  from 
this  step.  First  is  datci],  which  represents  the  unpolarized  portion  of  light  from  the  object  being  observed.  The  last 
three  data  sets  correspond  to  the  polarized  light  reflecting  off  the  object  in  observation.  Rectangle  five  represents  the 
process  for  adding  Poisson  noise  to  the  blurred  data  sets.  It  is  most  important  here  to  understand  that  dataj  is  used  to 
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calculate  the  three  unpolarized  noise  realizations  along  with  the  first  polarized  realization.  Finally,  the  simulated 
measured  data  is  calculated  by  adding  each  respective  unpolarized  component  with  its  polarized  partner.  The  three 
resulting  data  sets  and  are  shown  as  an  example  in  Figure  3.5. 


Figure  3.4:  Logical  flow  of  data  degradation  simulation  process.  Rhombuses  represent  creation  or  input  of  data  to 
simulator.  Rectangles  represent  calculations  made  in  the  simulation. 


polarized  blurred  data  w/  C2 


10  20  30 


polarized  blurred  data  w/  C3 
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Figure  3.5:  Examples  of  2  binary  stars  degraded  by  Poisson  modeled  noise  and  blurring  from  .  Left- 
Middle-  .  Right- 

Equations  (3.8-3.10)  are  the  average  value  of  and  respectively.  These  terms  can 


also  be  referred  to  as  intensities. 


and  .  The  intensity  measurements  look  very  similar  to  the 


measured  data,  but  have  no  noise  added.  An  example  of  an  intensity  image  is  shown  in  image  B  in  Figure  3.1. 
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Equation  (3.11)  shows  the  mathematical  relation  between  measured  data  and  intensity.  For  further  clarification  of 
what  the  other  variables  in  the  Equations  (3.8-3.10)  are  equivalent  to  refer  to  Table  3.1.  These  equations  are  a  simple 
consequence  of  Malus’  Law,  which  describes  the  intensity  of  the  light  passed  by  a  polarization  analyzer.  )  [4:332, 
13:1067].  The  discrete  convolutions  are  approximations  of  the  continuous  case.  These  equations  are  also  further 
explained  in  the  derivation  section  (3.2.2). 


(3.8) 


- 

(3.9) 

- 

(3.10) 

(3.11) 

Table  3.1:  Equivalencies  between  Equations  (3.8-3.10)  and  variables  described  in  Figure  3.4 


Equations  3.8-3.10 

Figure  3.4 

Notes 

data  set 

data  set 

Degree  of  Polarization 

True  intensity 

A.K.A  object 

data  set 

Fourier  Transform: 

Next,  an  intuitive  explanation  of  Equations  (3.8-3.10)  may  be  helpful  in  preparation  for  the  deconvolution 
model  (section  3.2),  since  it  describes  a  derivation  with  these  equations  as  a  starting  point.  Remember  that  since  the 
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optical  system  is  linear  and  space  invariant,  the  convolution  of  the  input  with  a  PSF  can  be  generalized  to  the 


superposition  integral  as  shown  below  (also  in  Equation  (3.1)),  but  now  the  polarization  analyzers  must  be 
accounted  for  in  the  convolution  integral.  The  double  sum  of  multiplied  by  the  system  PSF 

is  a  convolution  and  it  is  this  function  that  explains  how  objects  in  the  observation  plane  are  seen  in  the  detector 
plane.  Also,  notice  that  P,  Cj,  C2,  and  C?  are  just  scalars  on  the  range  from  (0,1).  The  next  thing  to  notice  is  that  the 
double  sum  is  broken  up  into  two  components.  The  first  component  represents  the  unpolarized  portion  of  light, 
while  the  second  is  the  polarized  portion.  These  equations  show  there  is  a  combination  of  unpolarized  light  along 
with  the  potential  for  polarized  light,  depending  on  the  object  in  observation.  It  is  shown  in  Malus’  Law  that 
approximately  half  of  the  unpolarized  light  will  be  transmitted  past  the  polarization  analyzer  [4:332]. 


3.2  Deconvolution  Model 

This  section  addresses  the  restoration  process  for  transformation  of  measured  data  to  an  estimated 

image  that  resembles  the  true  intensity.  The  improvement  in  image  sharpness  is  entirely  a  result  of  the  measured 
data  being  run  through  a  deconvolution  algorithm.  An  extensive  discussion  of  the  mathematical  derivation  regarding 
the  deconvolution  algorithm  takes  place  in  this  section  after  the  simulation  description  is  provided.  The 
deconvolution  algorithm  works  as  displayed  in  Figure  3.6. 

In  simulation  initial  estimates  for  seven  parameters  must  be  created  before  the  deconvolution  can  start. 
These  seven  parameters  are  shown  at  the  top  of  Figure  3.6.  These  estimates  constantly  improve  inside  the 
deconvolution  loop  until  a  declared  amount  of  iterations  are  exceeded.  Although  other  criteria  can  be  used  to  stop 
the  algorithm  they  would  require  user  input  as  does  a  maximum  iterations  criteria.  Therefore  in  both  cases  there  will 
be  no  algorithmic  autonomy.  This  concept  is  explained  later  in  Section  6.3.  Definitions  for  the  initial  estimates  are  as 
follows: 

1.  r0  =  Fried’s  Seeing  Parameter  (same  value  used  in  )  (refer  3.1.1,  and  [2:429] ) 

2.  =  Overall  System  Optical  Transfer  Function  Estimate  (refer  to  Equation  (3.4)) 
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3.  =  Object  Estimate  (arbitrarily  set  to  l)(refer  to  Table  3.1) 

4.  =  Degree  of  Polarization  Data  (arbitrarily  set  to  1/2)  (refer  to  3.1.3,  3.2.2-Equation  (3.32),  and 
3.2.3) 

5.  =  Polarization  Attenuation  Factors  (calculated  from  arbitrary  set  value  )  (refer  to  3.1.3  & 

Table  3.2) 

6.  =  Intensity  estimates  (calculated  from  exact  value  and  arbitrary  values  for 

and  )  (refer  to  Equations  (3.8-10)) 

7.  =  Conditional  Expected  Value  of  Complete  Data  Estimates  (calculated  from 

exact  and  measured  data  values,  and  arbitrary  values  for  and  )  (refer  to 

3.2.2  Step  4  [13:4],  [1]) 


Figure  3.6:  Deconvolution  algorithm  flow  chart:  The  parallelogram  consists  of  the  seven  necessary  initial 
estimates.  The  large  rectangle  consists  of  all  the  calculations  in  the  deconvolution  loop.  The  diamond  represents 
the  decision  to  terminate  to  the  deconvolution  loop. 

3.2.1  The  Update  Equations 


Update  equations  are  the  key  to  this  entire  research  effort.  The  results  of  these  equations  are  responsible  for 
the  restoration  process’  level  of  success.  The  deconvolution  loop  in  Figure  3.6  shows  the  parameters  that  are 
iteratively  updated.  Listed  below  in  Equations  (3.12  -  3.16)  are  the  corresponding  update  equations  (excluding 
intensity  and  conditional  expected  value  of  complete  data,  refer  to  Equations  (3.8  -  3.10)  and  (3.25-28)  respectively). 
Note  that  for  the  sake  of  brevity  that  all  plane  coordinates  are  not  included  for  the  the  degree  of  polarization, 

the  measured  data  and  the  complete  data  .  Following  this  listing  of 

equations  is  the  math  necessary  to  arrive  at  these  update  equations.  Section  3.2.2  will  briefly  explain  how  each  of 
these  update  equations  are  determined  and  what  they  are  solving  for. 
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(3.12) 


(3.13) 


(3.14) 


(3.15) 


(3.16) 


3.2.2  Partial  3  Channel  Polarimeter  Derivation 

A  good  portion  of  this  research  effort  lies  in  the  mathematical  derivation  for  the  three  channel  polarimeter 
restoration  algorithm.  This  derivation  relies  on  the  results  of  years  of  research,  as  discussed  in  Chapter  2.  This 
section  explains  the  use  of  MAP  estimation  techniques  in  combination  with  the  EM  process.  Essentially  this 
derivation  explains  the  mathematics  necessary  to  take  a  blurred  and  noisy  image  and  transform  it  into  a  sharper 
looking  image. 

The  six  steps  listed  below  will  briefly  summarize  the  derivation  for  the  three  channel  polarimeter 
restoration  algorithm,  while  the  full  derivation  is  available  in  Appendix  A.  The  first  three  steps  below  define  the 
mathematical  models  for  the  types  of  data  considered  in  this  research,  while  the  last  three  steps  are  the  signature 
steps  in  the  Expectation  Maximization  process. 
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Table  3.2:  Definitions  of  the  variables  for  the  3  Channel  Polarimeter  Derivation. 


Detector  Plane  coordinates: 

Observation  Plane  coordinates: 

Incomplete  Data  (measured): 

Complete  Data  Unpolarized  Component 
(mythical): 

Complete  Data  CH  1  Polarized 
Component  (mythical) 

Complete  Data  CH  2  Polarized 
Component  (mythical) 

Complete  Data  CH  3  Polarized 
Component  (mythical) 

True  Intensity  (object): 

Simulated  Image  (intensity): 

System  Point  Spread  Function: 

Pixels,  or  image  size: 

True  Degree  of  Polarization: 

Polarization  Attenuation  Factors: 

1/2 

Step  1:  Obtain  a  Statistical  Model  for  the  Incomplete  Data  (Measured  Data) 

Recall  Equations  (3.8-10)  from  section  3.1.3  and  the  double  sums  over  the  values  ,  which  are  the 
coordinates  in  the  observation  plane  (outer-space).  Notice  that  the  intensities  have  coordinates  in  the  plane, 
also  known  as  the  detector  plane  (camera  detector).  These  three  equations  mathematically  explain  the  occurrence  of 
a  convolution.  Also  recall  that  the  measured  data  (equivalent  to  incomplete  data)  is  .  It 

is  called  incomplete  data  because  it  is  not  a  direct  measurement  of  an  object  in  observation,  whereas  the  complete 
data  represents  a  direct  measurement  of  an  object  in  observation  (refer  step  2).  The  average  value  of  the  incomplete 
data  shown  below  and  in  Equations  (3.8-3.10)  builds  upon  Schulz’s  model  by  accounting  for  polarization  effects, 
thus  adding  the  true  degree  of  polarization  and  polarization  attenuation  factors.  [13:  1067] 
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Step  2:  Invent  a  Set  of  Complete  Data  (Mythical) 


The  complete  data  is  a  mathematical  construct  designed  to  be  a  direct  measurement  of  the 

image  without  the  convolution  operation,  therefore  it  represents  the  object  in  observation  before  its  light  passes 
through  the  atmosphere.  The  convolution  function  is  the  mathematical  reason  for  the  occurrence  of  blurring. 
Equations  (3.17-3.19)  are  broken  up  into  two  parts.  The  first  part  is  a  double  sum  of  the  unpolarized  component  of 
the  complete  data,  and  the  second  part  is  the  double  sum  of  the  polarized  component  of  the  first  channel  of  data.  The 
complete  data  are  chosen  to  be  independent  Poisson  distributed  random  variables,  which  means  that  the  incomplete 
data  will  also  be  a  Poisson  distributed  random  variable.  This  is  due  to  the  statistical  property  that  the  summation  of 
Poisson  random  variables  is  a  Poisson  random  variable.  [13:  1067]  Note  that  the  mean  of  unpolarized  complete  data 
are  mathematically  equivalent  so,  .  This  concept  is  applied  in  the  next  step. 


(3.17) 


(3.18) 


(3.19) 
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Step  3:  Select  a  Statistical  Model  for  the  Complete  Data 


The  main  goal  in  selecting  a  statistical  model  is  ensuring  statistical  consistency  between  the  incomplete  and 
complete  data.  As  stated  in  the  last  step  the  complete  data  in  this  derivation  was  chosen  to  be  modeled  as  Poisson 
data,  because  the  sum  of  Poisson  random  variables  is  a  Poisson  random  variable  [13:  1067,  EQN  15].  Equations 
(3.20-3.23)  show  the  mean  of  the  complete  data  sets. 


Where  in  unpolarized  case 


(3.20) 


(3.21) 


(3.22) 


(3.23) 


Step  4:  Formulate  the  Complete  Data  Log-Likelihood 

The  complete  data  log-likelihood  function  is  formulated  assuming  that  the  complete  data  are  statistically 
independent.  Equation  (3.24)  shows  the  log  of  the  joint  probability  for  all  of  the  complete  data.  Once  again,  this 
equation  resembles  Schulz’s  log-likelihood  function,  except  for  the  notational  differences,  additional  channels  of 
data,  and  polarized  data  (  ).  [13:  1067-Equation  17] 


3.24) 
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Notice  the  addition  of  polarization  prior  to  the  end  of  Equation  (3.24).  The  main  reason  for 

adding  this  is  to  constrain  the  values  of  polarization  to  be  between  0  and  1,  which  is  “prior”  known  information 
[7:310  &  350]. The  PDF  chosen  to  model  this  prior  is  a  Super  Gaussian,  .  This  PDF 

was  chosen  because  it  approximates  a  uniform  PDF  on  (0,1)  for  high  values  of  n.  A  uniform  PDF  has  to  be 
approximated  because  of  discontinuities  in  its  derivative  [6:23]. 

Step  5:  Expectation  Step 

The  conditional  expected  value  of  the  complete  data  log  likelihood  must  be  computed  in  order  to  formulate 
the  EM  algorithm.  The  conditional  expectations  of  each  piece  of  complete  data  are  conditioned  on  the  incomplete 
data  and  the  estimates  of  ,  ,  ,  .and  .  They  are  shown  below  in  Equations  (3.25- 

3.28)  and  computed  using  results  from  Shepp  and  Vardi  [14:115  &  119],  [13:1067,  EQN  26].  The  resulting 
conditional  expectations  will  differ  from  the  results  derived  by  Sheep  and  Vardi,  and  also  from  Schulz’s  results  for 
conditional  expectations  [14],  [13].  This  is  strictly  due  to  the  fact  that  this  research  is  using  three  channels  of 
polarimetric  data. 
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(3.25) 


(3.26) 


(3.27) 


(3.28) 


The  expected  value  of  the  complete  data  log  likelihood  is  added  to  the  natural  logarithm  of  the  prior  to 
yield  the  Q  function  shown  below: 


(3.29) 


Step  6:  Maximize  the  Expected  Value 

This  step  requires  the  most  algebraic  manipulation,  since  the  derivative  of  with  respect  to 
,  ,  and  must  be  calculated,  and  then  set  equal  to  zero.  For  the  sake  of 

brevity  only  the  results  of  each  derivative  will  be  shown,  along  with  the  four  conditional  expectations.  A  more 
detailed  derivation  is  shown  in  Appendix  A. 


(3.30) 


(3.31) 


(3.32) 
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The  most  important  operation  in  the  deconvolution  loop  is  the  estimation  of  .  This  is  accomplished 

by  solving  for  the  roots  in  Equation  (3.32).  Once  has  been  estimated  the  object,  intensities,  and  complete 

data  expectations  can  be  estimated  (refer  Figure  3.6). 

3.2.3  Point  Source  Pixel  Spacing 

One  last  concept  needs  explanation  before  the  simulation  results  can  be  discussed  and  this  is  how  the  ratio 
of  the  lens  diameter  determines  the  dual  point  source  pixel  spacing  in  the  detector  plane  (refer  to  Figure  3.7  and 
Equation  (3.33).  If  an  observer  is  viewing  two  stars  in  space  that  have  an  angular  separation  radians  through  a 
meter  diameter  lens,  then  the  resulting  spacing  in  the  detector  plane  will  be  meters.  For  example,  assume 
simulation  parameters  equal  to:  ,  resulting  in  .  If  the 

lens  diameter  changes  to  then  must  change  accordingly  to  .  Assuming  that  the  camera  pixel 

detector  size  =  ,  then  a  lens  would  have  spacing,  whereas  the  lens  would  have 

spacing.  To  ensure  consistency  in  all  simulations  this  fact  needs  to  be  considered,  because  the  angular 
separation  of  the  stars  will  not  change,  but  the  simulation  parameter  will  change  to  or  ,  resulting 

in  different  values  for  pixel  spacing.  These  two  lens  sizes  were  chosen  because  they  are  the  actual  sizes  of  telescope 
lenses  at  the  Maui  Space  Surveillance  Site  (MSSS)  [16].  Use  of  these  lens  sizes  ensures  more  meaningful  results  for 
Maui  Air  Force  Research  Labs,  the  sponsor  of  this  research. 


Figure  3.7:  is  the  angular  separation  of  2  point  sources  in  space,  is  the  diameter  of  the  lens  of  a  telescope  and 
is  the  distance  or  pixel  separation  in  the  detector  plane.  The  angular  separation  does  not  change,  but  a 
change  in  will  determine  a  change  in 


(3.33) 
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Another  topic  related  to  the  point  source  pixel  spacing  is  proper  sampling  of  data.  To  assure 


proper  Nyquist  sampling  criteria  consider  Equations  (3.34)  and  (3.35)  [12:57].  Assume  the  following  parameters  for 
a  given  scenario:  Also,  assume  that  focal  length  is 

adjustable,  which  in  the  case  of  MSSS  telescopes  is  true.  Based  off  of  Equation  (3.35)  these  parameters  determine 
to  ensure  proper  sampling  in  space.  If  ,  then  - ,  giving  a  Nyquist  sampling 

frequency  of  - . 


(3.34) 

(3.35) 
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4.  Simulation  Results 


Simulation  results  are  important  for  learning  the  ideal  performance  of  the  three  channel  polarimeter.  It  is 
possible  to  accurately  measure  performance  because  the  final  object  estimate  can  always  be  compared  to  the  true 
intensity.  This  comparison  allows  for  a  total  squared  error  (TSE)  comparison  and  most  importantly  a  resolution 
comparison.  Resolution  comparisons  are  important,  because  the  goal  of  this  research  is  to  “see”  whether  the  three 
channel  approach  can  improve  upon  the  two  channel  approach. 

This  section  analyzes  the  simulation  results  from  two  types  of  simulations.  The  first  type  of  simulation  is  a 
dual  binary  star  simulation.  In  space  imaging  a  dual  point  source  image  could  be  the  result  of  viewing  two  stars 
close  to  each  other  in  an  angular  separation  sense.  The  second  type  of  simulation  is  executed  with  an  actual  satellite 
image. 


4.1  Dual  Point  Source 

The  dual  point  source  simulation  results  are  shown  in  Figures  4.1  and  4.2  below.  Both  figures  are 
composed  of  six  images.  The  top  left  image  in  Figures  4.1  and  4.2  are  the  true  intensity.  It  represents  two  point 
sources  with  pixel  spacing  equal  to  nine  pixels  for  the  lens  simulation  and  twelve  pixels  for  the 

case.  The  degree  of  polarization  data  is  not  shown,  but  in  all  dual  point  source  cases  the  left  point  was  fully 
polarized  and  the  right  point  was  half  polarized.  The  top  right  image  in  both  figures  is  the  final  object  estimate  for 
the  two  channel  case.  All  other  images  in  both  figures  are  final  estimates  for  the  three  channels  at  different 
values.  The  corresponding  TSEs  are  displayed  in  a  table  below  each  representative  figure.  The  TSE  was  calculated 
using  Equation  (4.1)  below.  The  comparisons  between  the  two  and  three  channel  algorithms  are  consistent  in  that 
the  total  intensity  estimated  in  each  image  is  conserved.  This  ensured  that  the  total  number  of  photons  in  the  true 
intensity  and  final  object  estimate  are  the  same. 


(4.1) 

The  figure  captions  also  define  some  other  very  important  simulation  parameters.  First  is  the  lens  size  in 
centimeters.  Second  is  the  seeing  parameter  r0  in  centimeters.  The  third  value  requires  a  little  explanation,  since  it 
differs  for  the  three  and  two  channel  results.  To  ensure  consistency  in  comparisons  the  signal  strengths  must  be 
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relatively  equal.  This  means  that  if  the  true  intensity’s  two  point  sources  are  represented  by  200,000  photon  counts 
in  each  white  pixel  in  the  two  channel  polarimeter,  the  three  channel  must  have  2/3  the  amount  of  photons,  because 
it  splits  the  light  signal  into  three  channels,  versus  two  channels.  The  fourth  value  is  the  point  source  pixel  spacing. 
Fifth  is  the  number  of  iterations  used  in  the  deconvolution.  The  sixth  number  is  the  angle  of  polarization  in  the 
observation  plane.  This  angle  directly  impacts  the  three  polarization  attenuation  factors  (refer  to  the  first  calculation 
rectangle  in  Figure  3.4  to  see  mathematical  correlation).  Lastly,  the  numbers  along  the  x-axis  and  y-axis  show  the 
size  of  the  image.  Simulations  for  the  lens  used  an  image  size  of  32x32  pixels.  Simulations  for  the 

lens  used  46x46  pixel  images. 


4.1.1:  160  cm  Lens  Results 


Fig.  4.1:  Results  for:  Lens  =  160  cm;  ro  =  15  cm;  Signal  Strength  =  134k/200k  photoelectrons;  Point  Source  Spacing 
=  12  pixels;  Iterations  =  6000.  Note:  all  Final  Estimated  Objects  are  for  3  CH  results,  except  top  right  image.  The 
TSE  for  the  2  CH  is  nearly  double  the  3  CH  estimates  and  is  shown  in  Table  4.1.  The  visual  reproduction  of  the  true 
intensity  is  close  to  the  same  for  all  deconvolutions. 

Table  4.1:  TSE  values  for  160  cm  lens  results 


TSE  Value 

2  CH 

7.32E10 

3  CH  (a=123) 

3.31E10 

3  CH  (a=199) 

3.25E1U 

3  CH  (<x=300) 

3.26E10 
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4.1.2: 


120  cm  Lens  Results 


Fig.  4.2:  Results  for:  Lens  =  120  cm;  ro  =  15  cm;  Signal  Strength  =  134k/200k  photoelectrons;  Point  Source  Spacing 
=  9  pixels;  Iterations  =  6000.  Note:  all  Final  Estimated  Objects  are  for  3  CH  results  at  differing  a  values,  except  top 
right  image.  The  TSE  for  the  2  CH  is  nearly  double  the  3  CH  estimates  and  is  shown  in  Table  4.2.  Visually  the  3  CH 
is  better. 


Table  4.2:  TSE  values  for  120  cm  lens  results 


TSE  Value 

2  CH 

7.13E10 

3  CH  (<x=123) 

3.12E10 

3  CH  (<x=199) 

3.13E1U 

3  CH  (<x=300) 

3.15E10 

Visually,  it  is  difficult  to  discriminate  if  the  three  channel  algorithm  provides  a  closer  reproduction  of  the 
true  intensity,  which  effectively  demonstrates  better  resolution.  The  TSE  values  in  Table  4.2  provide  a  numerical 
metric  for  comparing  performance  showing  the  three  channel  algorithm’s  superiority. 
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Estimated  Object  3CH  Estimated  Object  2CH 


Fig.  4.3:  Results  for:  Lens  =120  cm;  rO  =  15  cm;  Signal  Strength  =  134k/200k  photoelectrons;  Point  Source 
Spacing  =  9  pixels;  Iterations  =  448k/500k.  Note:  the  color  inversion  from  other  images  and  that  the  3CH 
performs  better  because  both  of  the  darkest  pixels  are  exactly  where  they  are  supposed  to  be,  which  is  pixel  (16, 
16)  and  (25, 16).  The  TSE  values  are  1.22E10  for  the  3  CH  and  1.62E10  for  the  2  CH. 

The  best  way  to  actually  determine  which  algorithm  works  better  is  to  analyze  the  TSE  values.  In  all 

simulation  cases  at  6000  iterations  the  TSE  values  for  the  two  channel  are  more  than  double  that  of  the  three  channel 

results.  Figure  4.3  shows  an  exception  to  this  statement.  In  this  figure  it  is  evident  that  the  three  channel  algorithm 

outperforms  in  visual  reconstruction  of  the  true  intensity,  because  the  darkest  spots  in  the  image  are  actually  located 

in  correct  pixel  locations  with  respect  to  the  true  intensity.  These  pixel  locations  are  and  ,  hence  the 

separation  for  the  lens.  It  is  also  evident  that  the  TSE  value  is  lower,  but  only  by  25%  this  time. 

The  TSE  value  is  important  in  this  research  because  it  provides  quantification  for  both  algorithms’  performances. 

4.2  Satellite  Image 

The  satellite  image  simulation  provides  a  better  sense  of  both  algorithms’  performance  on  a  more  complex 
image.  There  are  three  given  sets  of  data  needed  for  the  satellite  image  simulation.  Only  the  first  two  sets  of  given 
data  are  shown  in  the  left  and  middle  of  Figure  4.4,  while  the  third  data  set  (angles  of  polarization  data)  is  not  shown 
because  it  looks  identical  to  the  polarized  data  image.  The  left  image  is  the  true  intensity,  which  is  the  ideal  satellite 
image.  The  second  image  is  the  polarization  data,  which  represents  the  parts  of  the  satellite  reflecting  polarized  light. 
The  right  image  is  the  measured  data,  which  is  a  blurred  version  of  the  truth  data  due  to  the  simulated  system  PSF. 
Figure  4.5  shows  the  images  produced  by  the  two  and  three  channel  algorithms. 
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4.2.1: 


160  cm  Lens  Results 


Truth  Data  Polarized  Data  Measured  Data 


50  100  150  200  250  50  100  150  200  250  50  100  150  200  250 


V 


\ 


Figure  4.4:  Left:  True  intensity.  Middle:  Polarized  data.  Right:  Measured  data  (160  cm  lens,  r0=22  cm) 


Figure  4.5:  Left:  True  intensity.  Middle:  3  channel  estimate  with  4000  iterations.  Right:  2  channel  estimate  with 
4000  iterations,  ro  =  22  cm.  The  TSE  value  of  the  3CH  is  2.58E10,  which  is  approximately  60%  the  value  of  the  2  CH 
TSE  at  4.2  5E10.  This  is  consistent  with  the  dual  point  source  simulations.  The  visual  results  of  the  true  intensity 
reconstruction  are  near  identical. 


4.2.2:  120  cm  Lens  Results 


Truth  Data 


Polarized  Data 


Measured  Data 
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Figure  4.6:  Left:  True  intensity.  Middle:  Polarized  data.  Right:  Measured  data  (120  cm  lens  blur,  ro  =  15  cm,  that  is 
why  the  measured  data  appears  more  blurry  than  the  measured  data  for  the  160  cm,  ro  =  22  cm.) 
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Figure  4.7:  Left:  True  intensity.  Middle:  3  CH  estimate  with  4000  iterations.  Right:  2  CH  estimate  with  4000 
iterations,  ro  =  15  cm.  The  TSE  value  of  the  3CH  is  2.94E10,  which  is  approximately  60%  the  value  of  the  2  CH  TSE 
at  4.79E10.  This  is  consistent  with  the  dual  point  source  simulations.  The  visual  results  of  the  true  intensity 
reconstruction  are  near  identical. 

When  viewing  the  results  for  the  satellite  deconvolution  in  Figure  4.7  there  is  no  visible  difference  between 
the  results  of  both  algorithms.  Even  if  larger  sized  images  were  provided  there  would  be  no  visible  difference.  Even 
though  the  three  channel  polarimeter  will  provide  slightly  higher  resolution  in  the  deconvolution,  it  is  not  significant 
enough  to  make  a  visible  difference  in  the  satellite  image  deconvolutions.  In  this  case  the  TSE  for  the  three  channel 
algorithm  was  not  quite  half  the  value  of  the  two  channel  algorithm,  but  closer  to  60%  of  the  two  channel  value. 


Remember  that  the  dual  point  source  simulations  for  a  lens  used  a  nine  pixel  separation,  and 

twelve  for  the  lens.  Figure  4.8  shows  how  poorly  both  algorithms  perform  when  the  pixel  separation  is 

significantly  decreased.  Now,  consider  the  satellite  image,  where  it  would  be  necessary  to  be  able  to  deconvolve 
point  sources  with  separations  as  little  as  two  or  three  pixels  to  get  a  much  better  estimate  of  the  true  intensity.  If 
both  algorithms  can  barely  deconvolve  at  nine  pixels  of  separation,  then  how  could  they  ever  deconvolve  details  of 
the  satellite  image,  which  have  details  with  two  or  three  pixels  worth  of  separation.  In  conclusion  the  algorithms  do 
have  slight  differences  in  performance,  but  not  enough  difference  to  visually  distinguish  between  both 
deconvolutions.  The  only  way  to  distinguish  between  the  three  and  two  channel  results  is  to  compare  TSEs  and 
binary  star  results.  The  TSE  for  the  two  channel  is  still  close  to  double  the  value  of  the  three  channel  result.  The 
error  is  evenly  distributed  through  the  illuminated  portion  of  the  images. 
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Figure  4.8:  Comparison  of  the  3  and  2  channel  algorithm  at  6  pixel  separation  shows  poor  performance  by  both. 
The  top  images  are  the  results  for  a  120  cm  lens  simulation  with  ro=15  cm.  The  bottom  images  are  the  results  for 
a  160  cm  lens  simulation  with  ro=15cm. 
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5.  Laboratory  Imagery  Results 


5.1  Laboratory  Setup 

Analysis  of  laboratory  imagery  with  the  three  channel  algorithm  takes  this  research  one  step  further  in 
proving  its  success  and  efficacy.  However,  the  given  imagery  and  results  in  this  section  can  only  model  actual  data 
gathered  from  a  telescope.  The  laboratory  images  were  gathered  as  shown  in  Figure  5.1.  The  light  source  for  this 
experiment  was  a  light  emitting  diode  (LED).  Light  passes  through  the  multi-sized  triple  slit  pattern  to  the  lens, 
where  it  is  then  focused  through  a  polarization  analyzer  to  the  CCD  camera. 


Figure  5.1:  Laboratory  Experiment  Setup:  light  source  is  an  LED.  Light  passes  through  the  multi-sized  triple  slit 
pattern  to  the  lens,  where  it  is  focused  down  to  the  CCD  camera  through  a  polarization  analyzer. 


5.2  Laboratory  Image  Blur  and  Noise 


The  blurring  that  would  be  caused  by  atmospheric  turbulence  in  astronomical  images  is  accounted  for  in 
the  lab  images  by  defocusing  the  camera,  so  the  laboratory  images  appear  blurred  (refer  to  Figure  5.2).  This  is  a 
reasonable  way  to  simulate  the  atmospheric  turbulence  assuming  the  severity  of  defocus  is  actually  comparable  to 
the  severity  of  blurring  caused  by  the  atmosphere.  One  way  to  assure  a  reasonable  amount  of  blurring  due  to  defocus 
is  to  compare  the  OTF  of  the  lab  data  to  the  simulated  OTF  as  explained  in  section  3.1.1.  This  is  a  reasonable 
comparison  because  the  OTF  can  be  thought  of  as  a  low  pass  filter.  If  the  filter’s  (OTFs)  can  be  matched  to  each 
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other  by  matching  cutoff  frequencies  and  overall  shape  of  the  curve,  then  it  is  reasonable  to  assume  that  they  will 
have  the  same  effect  on  any  given  truth  data. 


50 
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50  100  150  200  50  100  150  200 

Figure  5.2:  Left:  in  focus  image  (True  Intensity)  of  triple  slit  taken  from  CCD  camera.  Right:  out  of  focus  image  of 
triple  slit  pattern  to  simulate  blur  in  an  image  caused  by  atmospheric  turbulence.  Note:  the  polarization  analyzer 
was  removed  for  this  picture. 

The  lab  data  PSF  was  created  from  the  standard  Richardson  Lucy  deconvolution  algorithm  given  in 
Equation  (5.1)  [11],  [13:1066].  However,  instead  of  estimating  new  values  of  the  object  using  old  estimates  of  the 
object  (  )  (where  is  the  Fourier  transform  of  the  measured  data,  is  the  Fourier  transform  of 

intensity,  and  is  the  complex  conjugate  of  the  OTF),  now  the  PSF  is  going  to  be  estimated  in  a  similar  way  as  the 
object  was.  This  is  shown  in  Equation  (5.2),  where  the  symbol  represents  the  correlation  function. 


(5.1) 

(5.2) 
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Figure  5.3:  A  comparison  of  Lab  generated  OTF  vs.  Simulated  data  OTF  (  ).  The  plots  are  nearly 

identical  if  the  side  lobes  of  the  top  plot  are  ignored,  which  means  that  these  OTF’s  are  a  match  leading  to  the 
same  effect  on  any  given  truth  data.  Cutoff  frequencies  for  both  plots  are  approximately  -10k  and  10k  nr1. 

The  OTF  of  the  lab  data  is  the  Fourier  transform  of  after  the  algorithm  in  Equation  (5.2)  has 
converged.  A  plot  of  this  OTF  versus  a  plot  of  from  section  3.1.1  is  shown  in  Figure  5.3.  was 

generated  using  a  lens  diameter  of  and  an  of  ,  which  are  the  same  parameters  used  in  sections  4. 1 .2 

and  4.2.2.  The  bell  shapes  of  the  plots  are  a  very  good  match  to  one  another  and  they  almost  completely  attenuate  all 
spatial  frequencies  before  and  after  the  values  of  and  respectively.  This  is  excluding  the  fact 

that  the  lab  data  OTF  has  side  lobes,  which  is  something  the  ideal  simulated  data  will  never  deal  with.  The  key  to 
assuring  that  these  two  plots  are  reasonably  comparable  is  to  ensure  that  the  truth  data  image  in  simulation 
is  the  same  as  the  images  measured  for  the  lab  data 

The  dominant  source  of  noise  is  still  be  assumed  to  be  photon  count  noise  as  it  was  for  the  simulations, 
because  the  thermal  noise  is  insignificant  when  compared  to  the  photon  count  noise.  Since  the  PSF  was  calculated 
using  a  Richardson  Lucy  algorithm,  it  must  be  assumed  that  the  noise  is  modeled  with  Poisson  statistics,  because  the 
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derivation  to  arrive  at  Equation  (5.1)  requires  a  Poisson  model  for  the  probability  of  a  photon  arrival  (refer  to 
Equation  (3.6))  (Derivation  Available  in  Appendix  B). 

5.3  Laboratory  Results 

The  results  in  Figure  5.4  show  an  overwhelmingly  impressive  performance  by  the  three  channel  compared 
to  the  two  channel  algorithm.  There  are  three  different  results  shown  for  the  two  channel  algorithm  at  the  bottom  of 
the  figure.  Each  two  channel  result  is  based  off  of  a  deconvolution  with  the  unpolarized  channel  of  measured  data 
along  with  either  the  0°,45°,  or  90°  data.  The  two  channel  algorithm  seemed  to  have  some  issues  in  noise 
amplification,  which  explains  the  white  spots  all  over  the  images.  The  suspected  reason  for  the  three  channel 
algorithms’  success  is  mostly  due  to  the  polarization  diversity  provided  by  three  channels  worth  of  polarized  data, 
compared  to  the  two  channels’  need  for  an  unpolarized  and  arbitrarily  polarized  channel  of  data.  This  concept  and 
limitation  is  further  explained  in  section  6.4. 
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Figure  5.4:  Lab  results  at  1000  iterations:  Top  left:  in  focus  image  (True  Intensity)  of  triple  slit  taken  from  CCD 
camera.  Top  middle:  Unpolarized  measured  data.  Top  Right:  3CH  deconvolution  estimate.  Bottom  Left:  2CH 
deconvolution  estimate  with  0°  data  and  unpolarized  data.  Bottom  middle:  2CH  estimate  with  45°  data.  Bottom 
right:  2CH  estimate  with  90°  data. 
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6.  Limitations 


6.1:  ratio 

The  most  detrimental  limitation  in  this  research  is  the  reality  that  the  ratio  will  always  determine  the  severity 
of  blurring.  This  fact  is  inescapable.  Recall  Equation  (3.2),  where  the  LEOTF  is  described,  and  Equation  (3.39), 
where  the  maximum  spatial  frequency  is  given.  Remember  the  in  the  exponent  of  Equation  (3.2)  and 

notice  how  the  lens  diameter  in  the  equation  below  can  be  substituted  in  for  the  numerator.  This  proves  that  the 
ratio  of  lens  diameter  to  the  seeing  parameter  determines  level  of  blur  in  an  image. 


6.2:  No  Blind  Deconvolution: 

A  second  limitation  is  that  the  three  channel  algorithm  is  not  a  capable  of  performing  blind  deconvolution. 
It  was  originally  in  the  scope  of  this  research  effort  to  include  blind  deconvolution  capability  in  the  three  channel 
algorithm.  In  simulation  the  estimated  equaled  the  true  value.  Therefore,  the  system  OTF  was  always  predicted 
perfectly.  For  deconvolution  of  the  laboratory  data  the  PSF  was  calculated  by  use  of  the  Richardson  Lucy 
deconvolution  algorithm.  The  impact  of  not  being  able  to  do  blind  deconvolution  is  not  too  severe,  since  in  the  case 
of  advanced  laboratories  the  seeing  parameter  can  be  measured  separately. 

6.3:  Lack  of  Algorithmic  Autonomy  (no  stopping  criteria) 

A  potential  source  of  contention  in  this  research  is  the  lack  of  a  stopping  criterion  for  the  deconvolution 
algorithm.  The  previous  two  channel  research  by  Capt  James  included  the  use  of  a  stopping  criterion  [6].  The 
criterion  was  based  off  of  a  total  squared  error  comparison  between  the  incomplete  data  and  the  estimated  intensity 
per  iteration  versus  the  variance  of  the  corresponding  incomplete  data  [6:29].  The  theory  behind  the  stopping 
criterion  is  explained  in  his  research  and  the  comparison  equation  is  given  below  as: 


(6.1) 
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The  main  purposes  of  the  criterion  were  to  add  autonomy  to  the  deconvolution  algorithm,  and  most  of  all  to 
avert  over-iteration  of  the  estimated  object.  The  issue  of  over-iteration  is  of  no  consequence  in  this  research  and  is 
proven  below,  while  the  lack  of  autonomy  is  a  potential  limitation,  but  not  of  severe  consequence.  It  is  shown  by 
simulation  (Figure  6.2)  that  the  total  squared  error  of  the  final  object  estimate  versus  the  true  intensity  continuously 
trends  downward,  and  that  the  likelihood  of  the  object  being  estimated  as  the  true  intensity  always  increases  (Figure 
6.1).  Theoretically  the  likelihood  should  always  increase.  There  are  slight  spikes  in  the  TSE  plot  showing  that  the 
TSE  does  not  actually  always  go  down,  even  though  it  trends  downward.  This  does  not  mean  that  the  algorithm  does 
not  work  correctly,  but  is  more  likely  a  numerical  issue  in  MATLAB  (simulation  program).  Numerical  issues  like 
this  can  be  attributed  to  round  off  errors  in  floating  point  numbers,  or  some  of  the  logic  added  to  the  deconvolution 
algorithm  to  avoid  dividing  by  zero. 


X  105 


Figure  6.1:  Likelihood  of  3CH  estimated  object  converging  towards  true  intensity  at  300k  iterations.  This 
simulation  was  run  with  120cm  diameter  lens  and  an  ro  =  15cm. 
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Figure  6.2:  Total  Squared  Error  of  3CH  estimated  object  versus  true  intensity  at  448k  iterations.  This  simulation 
was  run  with  120cm  diameter  lens  and  an  ro  =  15cm. 

6.4:  Diversity  of  Multi-Channel  Data  due  to  Polarization 

The  more  diversity  there  is  in  a  set  of  data,  whether  it  is  two  or  three  channels  worth  of  data  will  give  a 
better  deconvolution  with  both  the  two  and  three  channel  algorithms.  In  the  case  of  the  laboratory  data  the  three 
channel  significantly  outperforms  the  two  channel  because  of  this  fact.  Figure  6.3  below  shows  that  the  45°  data  is 
significantly  affected  by  the  setting  of  the  polarization  analyzer.  This  channel  of  data  is  completely  missing  the 
lowest  bar  worth  of  light,  because  it  was  blocked  by  the  polarization  analyzer  at  a  45°  setting.  The  two  channel 
algorithm  did  not  perform  nearly  as  well,  because  it  required  an  unpolarized  channel  of  data  along  with  a  polarized 
channel  (0°  used  in  comparison  in  Chapter  5),  where  there  wasn’t  enough  diversity  to  allow  for  a  better 
deconvolution.  If  there  were  less  diversity  in  the  three  sets  of  data,  then  the  potential  for  success  in  deconvolution  is 
significantly  hampered.  This  is  one  significant  downfall  for  multi-channel  polarimeter  based  deconvolution 
algorithms. 


Figure  6.3:  3  channels  worth  of  laboratory  data.  Shows  that  the  polarization  analyzer  set  at  45°  impacts  the  triple 
bar  image  by  blocking  light  from  the  3rd  bar  (or  lowest  bar). 
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7.  Conclusions 


Atmospheric  turbulence  will  always  have  a  degenerative  effect  on  imagery,  which  is  dependent  on  the 
earth’s  atmosphere  as  a  medium  for  light  to  travel.  Polarimetry  in  conjunction  with  deconvolution  has  proven  to  be 
an  effective  method  for  mitigating  these  negative  atmospheric  effects.  The  reason  polarimetry  is  so  effective  is  that 
it  creates  the  potential  for  diversity  in  imagery,  which  is  advantageous  for  the  three  channel  algorithm.  More 
importantly,  it  shows  that  the  three  channel  algorithm  always  outperforms  the  two  channel,  even  though  in  some 
cases  resulting  imagery  may  appear  to  be  the  same.  Simulations  showed  that  TSE  estimates  were  always 
significantly  better,  while  visual  comparisons  for  the  three  channel  were  always  better  or  the  same  as  the  two 
channel.  Conversely,  in  the  laboratory  comparison  the  three  channel  provided  a  far  superior  reconstruction  of  the 
true  image,  while  the  two  channel  had  apparent  issues  with  noise  amplification. 

Though  there  are  some  limitations  with  the  three  channel  algorithm  as  described  in  Chapter  6,  there  is 
definite  potential  for  future  research  to  be  done  in  continuation  with  the  developed  algorithm.  This  research  effort 
only  consists  of  results  from  a  laboratory  and  simulation  environment,  so  the  obvious  next  step  is  to  use  the  three 
channel  algorithm  with  actual  telescope  imagery.  There  also  needs  to  be  more  research  done  in  the  area  of 
developing  effective  stopping  criteria  for  the  deconvolution  algorithm,  since  the  only  stopping  criteria  used  is  a 
declaration  of  maximum  iterations.  Once  this  is  successfully  completed  the  algorithm  can  start  to  be  used  in  an 
operational  environment  to  actually  assist  the  DoD  in  space  situational  awareness.  There  is  also  potential  for  the 
three  channel  algorithm  to  be  combined  with  concepts  of  other  research,  such  as  described  in  Milo  Hyde’s  IEEE 
article,  “Material  Classification  of  an  Unknown  Object  Using  Turbulence-Degraded  Polarimetric  Imagery.”  If  this 
research  is  continued  and  completed  it  may  provide  significant  capability  to  the  DoD,  and  as  a  result  better  national 
security  for  our  wonderful  United  States  of  America. 
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Appendix  A:  3  Channel  Polarimetric  Data  Deconvolution  Algorithm  Derivation 


Detector  Plane  coordinates: 

X,y 

Observation  Plane  coordinates: 

Z,W 

Incomplete  data  (measured): 

di(x,y),  d2(x,y),  d3(x,y) 

Complete  Data  mythical: 

True  unaberrated  image  (object): 

o(z,w) 

Simulated  image  (intensity): 

ii(x,y),  i2(x,y),  i3(x,y) 

Point  Spread  Function: 

h(x-z,  y-w) 

Pixels,  or  image  size: 

N 

Polarization  parameter: 

P(z,w) 

Transmission  of  3  polarizers: 

1/2 

Step  1:  Obtain  Statistical  Model  for  the  Measured  Data 


Note:  only  half  of  the  light  is  unpolarized,  since  the  polarizer  cuts  half  out 

(1) 

- 

(2) 

- 

(3) 

A -42 


(4) 


Step  2:  Invent  a  set  of  Mythical  Data  ( Complete  Data) 


(5) 

(6) 

Where  unpolarized  complete  data  is  the  same  so. 

(7) 

Step  3:  Select  a  Statistical  Model  for  the  Complete  Data 

Note:  data  chosen  to  be  Poisson  because  the  sum  of  Poisson  random  variables  is  a  Poisson  random  variable 


Where  P( z,w)=0  in  unpolarized  case 

(8) 

(9) 

(10) 
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(11) 


Step  4:  Formulate  the  Complete  Data  Log -Likelihood 


(12) 

(13) 

(14) 
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(15) 

(16) 

(17) 

(18) 

MAP  (Maximum  a  Posteriori)  Estimation:  MAP  Estimation  is  used  because  the  polarization  information  is  known  prior  to  maximization  calculations. 
MAP  Estimation  will  try  to  maximize  ,  while  Maximum  Likelihood  tries  to  maximize 


(19) 


can  be  disregarded  because  it  does  not  have  the  parameters:  P(z,w),  o(z,w),  or  C(z,w),  and  vanishes  when  the  derivative  is  taken  in  the 
Maximization  step. 


(20) 


The  higher  the  n  value  in  the  Super  Gaussian  PDF,  the  more  it  models  a  step  function.  In  reality  P(z,  w)  must  always  be  less  than  1.  This  PDF  is  chosen 
because  it  can  be  modeled  as  step  function  to  drop  at  one,  therefore  mathematically  making  it  very  unlikely  to  get  a  number  greater  than  1  for  P(z,w). 
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Step  5:  Expected  Value 


(21) 


Aside:  d  —  dl+d2,  where  t/is  true  measured  data,  dl=  one  instance  of  a  given  sum,  d2  =  all  other  instances  of  the  sum  c/is  a  sum  of  measured  data. 
Now  equate: 


Refer  to  EQNs  1-3  for  E[t/] 


(22) 


(23) 


(24) 

(25) 

(26) 

(27) 

(28) 
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(29) 

(30) 

(31) 

(32) 

(33) 

(34) 

(35) 

A -47 


(36) 

(37) 

Substitute  EQNs  29-32  into  38-41  respectively: 


(38) 

(39) 

(40) 

(41) 
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Notice  that  the  terms  dropped  off  in  EQNs  33-36.  This  happens  because  their  derivatives  will  be  equal  to  0  in  the  next  step,  since 

derivatives  will  be  with  respect  to  o(z,w),  not  o°,d(z,w). 


Derivative  of  natural  log  note: 


Rule  of  PSFs  note: 


will  always  double  sum  to  =  1,  because  the  law  of  the  conservation  of  energy  is  broken  if  the  double  sum  is  greater  than  1. 
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(47) 

(48) 

(49) 

(50) 

(51) 

(52) 
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(53) 


Algebra  to  simplify  o(z0,w0): 


Note:  h  =1  due  to  the  (refer  to  Rule  of  PSF’s  note  above  EQU  47) 

(54) 

-  - 

(55) 

(56) 

Now  differentiate  only  Q4  with  respect  to  C3(zo,w0)' 


(57) 

(58) 

Now  differentiate  Q  with  respect  to  P(zo,Wo): 


(59) 
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(60) 

(61) 

(62) 

Canceling  out  Ci  and  C2  terms: 


(63) 

(64) 

A -52 


Appendix  A:  3  Channel  Polarimetric  Data  Deconvolution  Algorithm  Derivation 


2011 


(68) 

(69) 

Last  Step:  Substitute  EQU  56  into  EQU  68 


-  -  -  _ 

(70) 

(71) 

Note:  Ignore  common  denominator,  because  the  main  idea  of  maximization  after  taking  the  derivative  of  all  the  Qterms  is  to  set  - equal  to  zero.  This 

happens  when  the  numerators  are  equal  to  zero. 

(72) 

(73) 

(74) 
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Appendix  B:  Derivation  for  Assumed  Poisson  Noise  in  Richardson  Lucy  Algorithm 


(Provided  by:  Dr.  Stephen  Cain,  ENG780-Statistical  Optics  Class  Notes) 


Deconvolution  with  Poisson  Noise 


N  N 
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i(z,w)dM  e~iiz'w) 


P[d(z,w)]  = 
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Deconvolution  with  Poisson  Noise 

L(o)  =  tt  d  (z,  w)  In  (i(z,  w))  -  i(z,  w) 

Z=1  w= 1 


dL(o) 


N  N 


■=zr 
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Deconvolution  with  Poisson  Noise 
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